{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "# Homework 1\n",
    "\n",
    "1. Write all code in the notebook.\n",
    "1. Write all text in the notebook. You can use MathJax to insert math or generic Markdown to insert figures (it's unlikely you'll need the latter). \n",
    "1. **Execute** the notebook and **save** the results.\n",
    "1. You can save the notebook as a PDF as well.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-01-22T19:57:47.188990Z",
     "start_time": "2019-01-22T19:57:46.107420Z"
    },
    "slideshow": {
     "slide_type": "-"
    }
   },
   "outputs": [],
   "source": [
    "%mavenRepo snapshots https://oss.sonatype.org/content/repositories/snapshots/\n",
    "\n",
    "%maven ai.djl:api:0.8.0\n",
    "%maven org.slf4j:slf4j-api:1.7.26\n",
    "%maven org.slf4j:slf4j-simple:1.7.26\n",
    "\n",
    "%maven ai.djl.mxnet:mxnet-engine:0.8.0\n",
    "%maven ai.djl.mxnet:mxnet-native-auto:1.7.0-backport"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import ai.djl.Device;\n",
    "import ai.djl.ndarray.NDManager;\n",
    "import ai.djl.ndarray.NDArray;\n",
    "import ai.djl.ndarray.LazyNDArray;\n",
    "import ai.djl.ndarray.index.NDIndex;\n",
    "import ai.djl.ndarray.types.Shape;\n",
    "\n",
    "// You have to use the MXNet engine for Q1 since\n",
    "// it is currently the only engine in DJL\n",
    "// to support Lazy NDArrays\n",
    "NDManager manager = NDManager.newBaseManager(Device.cpu(), \"MXNet\");"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%load ../../utils/StopWatch"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Speedtest for vectorization\n",
    "\n",
    "Your goal is to measure the speed of linear algebra operations for different levels of vectorization. You need to use the LazyNDArrays `waitToRead()` function on the output to ensure that the result is computed completely, since NDArray uses asynchronous computation. Please see https://javadoc.io/doc/ai.djl/api/latest/ai/djl/ndarray/LazyNDArray.html for details. \n",
    "\n",
    "Hint: The MXNet NDArray manager utilizes MxNDArrays which implement the LazyNDArray class. \n",
    "You must cast the return NDArray to LazyNDArray to be able to call\n",
    "the `waitToRead()` function.\n",
    "\n",
    "1. Construct two matrices $A$ and $B$ with Gaussian random entries of size $4096 \\times 4096$. \n",
    "1. Compute $C = A B$ using matrix-matrix operations and report the time. \n",
    "1. Compute $C = A B$, treating $A$ as a matrix but computing the result for each column of $B$ one at a time. Report the time.\n",
    "1. Compute $C = A B$, treating $A$ and $B$ as collections of vectors. Report the time.\n",
    "1. Bonus question - what changes if you execute this on a GPU?\n",
    "\n",
    "Note: The vector by vector calculation will likely take over 1.5 hours on a modern cpu."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "var A = manager.randomNormal(new Shape(4096, 4096));\n",
    "var B = manager.randomNormal(new Shape(4096, 4096));\n",
    "var sw = new StopWatch();\n",
    "var C = (LazyNDArray) A.dot(B);\n",
    "C.waitToRead();\n",
    "System.out.printf(\"Matrix by matrix: %f seconds\\n\", sw.stop());\n",
    "\n",
    "C = (LazyNDArray) manager.zeros(new Shape(4096, 4096));\n",
    "sw.start();\n",
    "for (int i = 0; i < 4096; i++) {\n",
    "    C.set(new NDIndex(\":\").addIndices(i),\n",
    "            A.dot(B.get(new NDIndex(\":\").addIndices(i))));\n",
    "}\n",
    "C.waitToRead();\n",
    "System.out.printf(\"Matrix by vector: %f seconds\\n\", sw.stop());\n",
    "\n",
    "// Note: this will likely take over 1.5 hrs on a standard cpu.\n",
    "C = (LazyNDArray) manager.zeros(new Shape(4096, 4096));\n",
    "sw.start();\n",
    "for (int i = 0; i < 4096; i++) {\n",
    "    for (int j = 0; j < 4096; j++) {\n",
    "        C.set(new NDIndex(i, j),\n",
    "                A.get(new NDIndex(i).addIndices(\":\"))\n",
    "                        .dot(B.get(new NDIndex(\":\").addIndices(j))));\n",
    "    }\n",
    "}\n",
    "C.waitToRead();\n",
    "System.out.printf(\"Vector by vector: %f seconds\\n\", sw.stop());"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Semidefinite Matrices\n",
    "\n",
    "Assume that $A \\in \\mathbb{R}^{m \\times n}$ is an arbitrary matrix and that $D \\in \\mathbb{R}^{n \\times n}$ is a diagonal matrix with nonnegative entries. \n",
    "\n",
    "1. Prove that $B = A D A^\\top$ is a positive semidefinite matrix. \n",
    "1. When would it be useful to work with $B$ and when is it better to use $A$ and $D$?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "1. $B$ is positive semidefinite if \\forall $x \\in \\mathbb{R}^m, x^TBx \\geq 0$. Then let $x \\in \\mathbb{R}^m$ and let $y = A^Tx$. Then\n",
    "$$x^TBx = x^TADA^T = (A^Tx)^TDA^Tx=y^TDy=\\sum_{1}^{n}d_iy_i^2 \\geq 0$$\n",
    "1. It would be more useful to work with $B$ if $m << n$ as matrix multiplication for instance for two arbitrary matrices $X$ and $Y$ with dimensions $a$ by $b$ and $b$ by $c$ is a runtime of $O(abc)$. \n",
    "If we multiply $B$ by a matrix $C$ that is $m$ by $k$, then $BC$ is computed in $O(m^2k).$ \n",
    "This matrix multiplication with $ADA^TC$ is $O(2mnk + n^2k)$. \n",
    "Thus if $m << n$, then it would be more efficient to use $B$ and if $n >> m$, then it would be better to use $A$ and $D$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. MXNet on GPUs\n",
    "\n",
    "1. Install GPU drivers (if needed)\n",
    "1. Install MXNet on a GPU instance\n",
    "1. Display `!nvidia-smi`\n",
    "1. Create a $2 \\times 2$ matrix on the GPU and print it. See https://d2l.djl.ai/chapter_deep-learning-computation/use-gpu.html for details."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "var managerGPU = NDManager.newBaseManager(Device.gpu());\n",
    "var x = managerGPU.ones(new Shape(2, 2));\n",
    "x;"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Memory efficient computation\n",
    "\n",
    "We want to compute $C \\leftarrow A \\cdot B + C$, where $A, B$ and $C$ are all matrices. Implement this in the most memory efficient manner. Pay attention to the following two things:\n",
    "\n",
    "1. Do not allocate new memory for the new value of $C$.\n",
    "1. Do not allocate new memory for intermediate results if possible."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "var A = manager.randomNormal(new Shape(100, 100));\n",
    "var B = manager.randomNormal(new Shape(100, 100));\n",
    "var C = manager.randomNormal(new Shape(100, 100));\n",
    "\n",
    "C.addi(A.dot(B));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. Broadcast Operations\n",
    "\n",
    "In order to perform polynomial fitting we want to compute a design matrix $A$ with \n",
    "\n",
    "$$A_{ij} = x_i^j$$\n",
    "\n",
    "Our goal is to implement this **without a single for loop** entirely using vectorization and broadcast. Here $1 \\leq j \\leq 20$ and $x = \\{-10, -9.9, \\ldots 10\\}$. Implement code that generates such a matrix."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "var x = manager.arange(-10, 10.1f, 0.1f).reshape(new Shape(201, 1));\n",
    "var j = manager.arange(1, 21, 1f).reshape(new Shape(1, 20));\n",
    "x.pow(j);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "celltoolbar": "Slideshow",
  "kernelspec": {
   "display_name": "Java",
   "language": "java",
   "name": "java"
  },
  "language_info": {
   "codemirror_mode": "java",
   "file_extension": ".jshell",
   "mimetype": "text/x-java-source",
   "name": "Java",
   "pygments_lexer": "java",
   "version": "11.0.5+10-LTS"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
